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ABSTRACT 


A summary  of  the  general  theory  of  nonlinear  differential  equations 
of  first  order  is  given,  with  the  aim  of  providing  practical  working  rules 
for  the  analysis  of  technical  problems,  without  pretense  to  rigor  and  com- 
pleteness. In  general,  only  equations  of  first  order  are  considered  here. 
After  a discussion  of  the  existing  conditions  and  the  analysis  of  singular 
points  and  of  the  integral  solutions  of  the  few  types  of  equations  which  can 
be  integrated  in  closed  form,  principal  analytical  and  graphical  procedures 
for  the  approximation  of  the  solutions  are  described. 
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INTRODUCTION 

The  analysis  of  the  d°c  response  of  single  energy  nonlinear  circuit* 
va»  presented  in  a previous  report.*  It  was  shown  that  rigorous  solutions  in 
closed  fora  could  be  obtained,  and  it  was  outlined  that  in  inis  case  sereral 
concepts  and  characteristics  derived  froa  linear  analysis  do  not  apply* 

It  is  nov  proposed  to  extend  the  analysis  to  the  a-c  response  of  the 
sane  circuits.  This,  however,  requires  the  application  of  special  techniques 
of  nonlinear  analysis  and  presents  difficulties  of  such  higher  order  than  the 
ones  previously  aet.  In  general,  no  rigorous  solutions  in  closed  fora  can  be 
obtained  and  one  has  to  sake  recourse  to  approximation  aethods.  Tor  this 
reason,  before  proceeding  to  the  analysis  of  the  «-c  response,  which  will  be 
presented  in  a future  report,  it  is  convenient  to  suisaarise  briefly  the  prin- 
cipal aethods  of  analysis  of  linear  and  nonlinear  differential  equations  which 
describe  the  behavior  of  single  energy  syrteas.  The  following  notes,  which  are 
—Inly  Halted  to  differential  equations  of  first  order,  are  only  intended  to 
give  working  rules  to  the  student  of  nonlinear  problems . The  reader  interes- 
ted in  rigorous  and  complete  discussion  of  the  subject  is  referred  to  the 
various  excellent  texts  on  differential  equations  published.** 


I.  General  Characteristics  of  Difference!  Equations  of  First  Order:  Exist- 

ence of  Solutions,  Regular  and  Singular  Points 


type 


A group  of  functions  represented  analytically  by  an  expression  of 


F (x,y,c)  • 0 


(1) 


where  c is  a variable  paraaeter,  is  said  to  constitute  a "family*.  The  charac- 
teristic relationship  of  the  functions  of  the  family  is  an  equation  obtained 
by  elimination  of  the  pars— ter  c between  (l)  and  its  derivative 


%r 

>x 


♦ 


*F  dv 

>y  <5 


0 


(2) 


Report  R-271-52,  PIB-210.  This  is  the  second  of  a series  of  reports 
on  the  analysis  of  nonlinear  circuits,  bassd  on  the  — terial  of  the 
Graduate  Course  "Nonlinear  Analysis*  offered  by  Dr.  B.  Saber  u t the 
Polytechnic  Institute  of  Brooklyn. 

See  for  instance:  X,  Picard,  Traits1  d* analyse,  Gauthier-Yillars 

(1905)*  H,  Poincare’,  Les  aethodes  nouvellee  de  la  aechanique  celeste, 

X*  Flsaaaricn,  Peris  (1908);  Z.  Goursat,  Differential  equations  (Sag. 
Trans.)  Boston  (1917);  B.  L.  Ince,  Ordinary  differential  equations, 
Dover  Fubl,  New  fork  (1927);  L.  Bieberbach,  Differentialgleichungen, 
Dover  Publ.  (1930);  G.  Seasons,  Bqpasioni  differensiali  nel  caapo  reals, 
Zanichelli,  Italy  (191*0). 
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ouch  relationship  of  type 

* (*,y,p)  “0  (3) 

where  p « dy/dx,  ia  called  a differential  equation  of  first  order.  Similarly, 
a family  of  functions  with  two  variable  independent  parameters  is  represented 
by  a differential  equation  of  second  order.  Therefore,  in  general,  the  order 
of  the  differential  equation  is  equal  to  the  number  of  arbitrary  parameters 
upon  which  the  family  of  functions  depend,  provided  that  these  parameters  are 
not  mutually  interdependent. 

The  reverse  process  of  finding  the  various  functions  of  the  family 
from  their  differential  equation  is  usually  very  difficult.  In  addition,  it 
is  clear  that  the  differential  equation  may  possess  other  integrals  besides 
those  represented  by  the  family  from  which  it  has  been  derived. 

For  example,  the  circles  in  a plane  farm  a three-dimensional  family 
of  equation 

x2  ♦ y2  ♦ 2Ax  ♦ 2By  ♦ C - 0 (U) 


Differentiating  three  times  and  eliminating  A,  B,  and  C,  one  has 

in  1 2 i 1 1 2 

y (l+y  ) - 3 y y - o 


This  differential  equation  is  satisfied  not  only  by  Eq.  (U),  but  also  by  the 
equation  of  any  straight  line  in  the  plane,  since 


• • 


- y 


■ * • 


- o 


is  a solution. 


In  the  following  we  shall  investigate  the  problem  of  the  existence 
of  solutions  in  a certain  domain  D(x,y)  (such  a domain  is  four-dimensional  if 
x and  y are  considered  complex  variables),  the  characteristics  of  various  types 
of  singular  points,  and  the  methods  of  solution  of  the  differential  Eq.  (3) 
with  closed  form  or  with  approximate  or  numerical  expressions. 


Regular  and  Singular  Points 

Given  a differential  equation  of  type  (3)  and  a domain  D(x,y)  one 
is  confronted  with  the  problems  of  investigating  a)  the  existence  of  solu- 
tions of  (3)  at  each  point  of  D,  b)  the  nature  of  such  solutions,  i.e.  whether 
or  not  they  all  belong  to  a family  of  type  (1),  and  c)  the  characteristics  of 
singular  points  of  D(x,y). 

Let  us  assume  that  Eq.  (3)  can  be  solved  for  p - dy/dx  and  that,  for 
any  initial  pair  (xo,y0)  in  D,  there  exists  one  and  only  one  root 
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P - f (x,y)  (5) 

which  reduces  to  Pq  when  r * Xq  ,y  - y0.  For  greater  generality  we  A all  con- 
sider x and  y as  complex  variables.  Cauchy  proved  the  following  existence 
theorem.  If,  within  a circle  |x-xcj<h,  f(x,y)  is  analytic  in  x and  in  y, 
then  Eq.  (5)  possesses  a unique  solution  y - y(x)  which  is  analytic  within 
the  circle  and  reduces  to  y0  when  x - Xo.  In  the  case  of  real  rariables, 
Lipschits  simplified  the  proof  of  Cauchy,  showing  that  the  existing  conditions 
reduce  to  a)  the  continuity  of  f (x,y)  within  the  rectangular  doaain 

|x-x0|<h  }y-y0Ub 

where  h:£  b/K  and  M is  the  upper  bound  of  ff(x,y)  | in  the  domain,  and  b)  the 
existence  of  a positive  number  k,  such  that 

|f(x,y’)  - f(x,y)|  < k|y*  - y| 

where  y and  y are  any  two  numbers  of  the  rectangular  domain. 

The  previous  existence  theorem  is  not  applicable  if  Sq.  (3)  possesses 
a multiple  p root  for  x - Xq,  y « y0.  In  this  case,  Eq.  (3)  is  equivalent  to 
a differential  equations  (where  m is  the  order  of  multiplicity  in  p)  and,  in 
ganaral,  possesses  a integrals  at  XqJo  called  singular  integrals.  Such  a 
situation  arises  at  points  on  the  envelope  of  the  family  of  integral  curves 
(1),  at  aultipla  points  of  any  integral  curve  (1)  where  two  branches  of  the 
same  curve  touch,  and  at  tac  points  where  two  no nconsecutive  curves  (1)  touch. 
The  totality  of  loci  for  which  at  least  two  values  of  p are  equal  is  obtained 
eliminating  p between  (3)  and  the  equation  d$/dp  - 0.  If  B(x,y)  • 0 is  the 
equation  so  obtained,  this,  in  ganaral,  will  not  satisfy  tha  given  differen- 
tial aquation  and,  for  this  reason,  will  not  be  one  of  its  integrals.  In  such 
casa,  R(x,y)  - 0 is  the  locus  of  the  multiple  points  or  the  locus  of  the  tac 
points  of  the  actual  integral  curves.  If  R(x,y)  - 0 satisfies  Eq.  (3),  then 
it  represents  the  envelope  of  the  family  of  integral  curves. 

To  illustrate  graphically  this  result,  one  can  consider  the  family 
of  isoclinals  obtained  from  Eq.  (3)  by  letting  p - const.  • X.  A general 
survey  of  the  integral  solutions  of  Eq.  ( 3)  is  obtained  by  sketching  such  n 
family  as  a function  of  the  parameter  k.  If  the  coordinate  plane  (x,y)  is 


Fig.  1 - Envelope  (E)  of  the  family  of  isoclinal  curves  and  locus  of  cusps 
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covered  with  enough  curves , it  is  possible  to  sketch  integral  curves  starting 
at  any  initial  point  (xo,y0)  and  proceeding  in  steps  from  (x,y)  to  x*dx,y*3jr 
in  the  direction  of  the  isoclinal  passing  through  (x,y).  In  Fig.  1 is  repre- 
sented the  envelope  S of  isoclinal  curves  and  indicated  that  it  is  the  locus 
of  cusps  of  the  integral  curves.  As  a matter  offset,  in  general,  the  iso- 
clinal curves  have  a slope  different  than  the  values  of  k for  which  they  are 
defined.  There  follows  that,  in  general,  integral  curves  cross  the  isoclinal 
curves.  However,  at  the  envelope  E,  which  also  has  a slope  different  than  p, 
the  integral  curves  cannot  cross  since  there  are  no  contiguous  isoclinal  curves 
on  the  other  side  of  the  envelope.  Consequently,  the  integral  curves  have  a 
cusp  or  a stop  point  on  the  envelope. 

The  equation  R(x,y)  - 0 can  also  include  a locus  of  double  points  of 
the  isoclinal  family.  The  case  is  indicated  in  Fig.  2 where  it  appears  that 
the  corresponding  curve  is  the  locus  of  tac  points  of  the  integral  curves.  As 


lig.  2 - Tac  locus  (F)  of  integral  curves,  and  locus  of  nodes  of  isoclinal 
curves. 


a natter  of  fact,  for  each  point  of  the  curve  F there  are  two  possible  direc- 
tions of  the  corresponding  integral  curves,  and  therefore  F is  the  locus  of 
tangency  where  2 nonconsecutive  curves  touch . 


In  order  to  find  the  locus  of  multiple  points  of  the  lsoclinals,  it 
is  necessary  and  sufficient  that  the  slope  of  the  isoclinals,  defined  by 

1*  * & - 0 
x )y  dx 


is  satisfied  by  two  different  values  of  dy/dx  if 


0* 

Ox 


0, 


0. 


Finally,  let  us  consider  (x,y)  pairs  for  which  Eq.  (3)  does  possess 
a unique  root  p,  but  the  existence  conditions  are  not  satisfied.  Such  pairs 
are  called  singular  points  of  the  equation  and,  in  general,  are  isolated. 


:Si 

'&  # 


i *■ 

1 


,4 

.J 
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They  play  an  important  part  in  characterising  the  behavior  of  the  various 
integral  curves  in  their  vicinity  because  they  ere  themselves  singular  points 
of  the  solutions,  i.e.,  points  at  which  the  solutions  are  discontinuous,  not 
unique,  or  not  existent. 

Let  us  assume  that  in  Eq.  (£),  f(x,y)  is  a rational  function  of  y. 


a . r(z,j)  . hill  (6) 

dx  Q(x,y) 

where  P(x,y),  Q(x,y)  are  polynomials  in  x,y  with  no  common  polynomial  divisor 


P(x,y)  - Pjx)  ♦ yP1(x)  ♦ y^Pg(x)  ♦ 


♦ A„(x) 


<Xx,y)  - QqCi)  ♦ yi^fx)  ♦ yz%(x)  * 


y\(x> 


The  singularities  of  (6)  are  discrete  and  may  be  separated  into  two  fundamental 
classes,  the  fixed  or  intrinsic  and  the  movable  or  parametric  singularities. 

The  first  ones  arise  at  points  x,  where  a)  any  of  the  coefficients  or  % has 
a singularity  which  cannot  lx  removed  by  multiplying  P(x,y)  and  Q(x,y)  times 
an  appropriate  function  of  xj  b)  Q(x,y)  is  identically  seroj  and  c)  tj(x^,y)«0, 
P(xi,y)  - 0 are  satisfied  simultaneously  by  a particular  value  of  y. 

The  Intrinsic  singularities  are  connected  with  only  some  of  the 
functions  F(x,y,c)  - 0 and  arise  in  correspondence  of  values  xi  where  such 
functions  have  a multiple  point.  For  this  reason,  they  depend  upon  the  value 
of  c,  i.e.,  upon  the  initial  conditions.  Parametric  singularities  are  found 
only  in  nonlinear  differential  equations;  they  can  coexist  with  intrinsic 
singularities. 

For  example,  the  sslution  of  the  differential  equation 

7--T2 

with  y(0)  - y0  is 

y - y0/U  - *y0) 

This  integral  has  a pole  at  x - l/jr0. 

It  may  be  shown  that  an  equation  of  the  first  order  and  second  de- 
gree of  the  Riccati  type 

& - p0(x)  ♦ yPi(x)  ♦ y2P2(x)  (7) 

dx 


&■ 

& 
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cannot  possess  multiple  point  singularities.  There  follows  that  linear  equa- 
tions of  first  order  and  nonlinear  Rlcoati  equations  possess  solutions  which 
are  rational  functions  of  the  Initial  value  y0>  i«e*> 

70  A(x)  * 3(x) 

y 

7q  C(x)  ♦ D(x) 

A Riccati  equation  can  be  converted  into  a linear  homogeneous  equation  of 
second  order  with  the  substitution 


y - - 


u P2(x) 


One  obtains  from  (7) 


2 

P2(x)  [.(X)  + P1(x)  P2(*)]  ^ «•  Po(x)  F*(x)u  - 0 

and  the  solution  of  this  differential  equation  is,  in  genera],  of  typo 

u - CjUj^x)  ♦ C^U). 

It  may  be  shown  that  the  movable  singular  points  of  an  equation  of 
type  (6)  can  only  be  poles  or  algebraic  critical  points. 


Classification  of  Intrinsic  Sine 


Points 


He  shall  now  examine  the  methods  for  the  determination  of  the  charac- 
teristics of  the  intrinsic  singular  points  of  Eq.  (6),  and,  in  particular, 
those  which  are  common  zeros  of  P(x,y)  and  Q(x,y),  but  are  not  stationary 
points  of  either  function.  If  (*l,yi)  is  one  such  pair,  expanding  P(x,y)  and 
<3(x,y)  in  its  vicinity  by  Taylor  series 

f(x»y)  -^(Xiyi)  +|^(x-x1)^(x1y1)  + (y-y^fyU^y^J  ♦ 


^ • o • e 

one  can  write  Eq.  (5)  as  follows: 

Painleve  - Equations  differentielles  ordinaires. 
t.  2,  V.3,  1910. 


Encycl.  des  Sciences  math. 


• x-  . 
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d£  _ &( x-x1)  ♦ bCy^)  ♦ P1(z,y)  } 

dx  cCr-^)  ♦ dCy-y^  ♦ Qj^y) 

In  Eq.  (8),  a - (*£)  , b -(i^)  , c - (2S)  , d - (|S)  , and  P^y), 

hx-A  *y  Vi  '6xxi^i  7 ¥1 

^(x,y)  are  polynomials  in  (x-x^),  (y-y1)  of  degree  not  less  than  two.  If  the 
differences  (x-x^) ,(y-y1)  are  made  to  approach  aero  along  some  arbitrary  curre 
in  (x,y),  P1(x,y)  and  Q1(x>y)  will  vanish  to  an  order  higher  than  the  first  so 
that,  in  a region  sufficiently  small  surrounding  (x^y^),  Eq.  (8)  may  be  written 
approximately 


£ - *<*-*!>  * 
dx  cU-x^  ♦ d(y-y^) 


(9) 


Eq.  (9)  is  homogeneous  in  (x-xx),  (y-yi)  and,  by  means  of  the  substitution 
y-yi  - s(x-x^),  may  be  reduced  to  one  whose  variables  are  separable.  Neglect- 
ing the  case  ad  - be-  0 for  which  equation  (9)  becomes  p - const.,  it  is 
found  that  the  characteristic  of  the  singularity  (x^^)  depends  upon  the  nature 
of  the  roots  *1,  *2  ot  the  algebraic  equation 

- (b*c)  % - ( ad- be)  - 0 
2 

i.e.,  upon  the  discriminant  l - (b-c)  ♦ U ad.  If  £-c0,  (x^y^)  is  a limit 
point  (focus)  of  the  integral  curves  which  arc  spiral-like}  if  A 2:0  but 
ad  - bc<0,  (xjjx)  is  an  actual  common  point  of  the  integral  curves  (node)} 
if  and  ad  - be  70,  (xjyi)  is  a saddle  point. 


Pig.  3 
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Letting 


-(Uc)  - Q - -(ad-bc)  - a1*2 


one  has 


A - P -UQ 


bo  that  the  previous  conditions  in  a P,  Q plane  correspond  to  the  regions  out- 
lined in  Pig.  3. 

It  is  possible  to  assign  an  index  to  the  various  types  of  singulari- 
ties. As  a natter  of  fact,  if  in  the  (xfy)  plane  we  consider  an  arbitrary 
closed  curve  C which  does  not  possess  any  multiple  points  and  surrounds  one 
and  only  one  singularity,  the  total  number  of  revolutions  K made  by  the  vector 
of  components  P(x,y),  Q(x,y)  in  a complete  circuitation  of  C is  +1  if  the 
singularity  is  a node  or  a focus,  -1  if  it  is  a saddle  point.  The  number  N 
is  called  index  of  the  singularity;  if  the  curve  C encloses  several  singulari- 
ties, the  value  of  If  is  the  sum  of  their  indexes.  This  nuaber  is  expressed 
mathematically  with  the  relation 


If  - — 


d tan' 


P(*,y)  2! 


FaQ  - QdP 

T7T 


As  a consequence  of  the  previous  statements  there  follows  that  if  the  curve  C 
is  an  actual  solution  of  Sq.  (9),  the  corresponding  value  of  N is  +1,  so  that 
the  sura  of  the  indexes  of  the  singularities  enclosed  must  add  up  to  +1. 

The  subject  of  singularities  of  a differential  equation  of  first 
order  will  be  discussed  further  in  a subsequent  report,  to  consider  cases  in 
which  x and  y are  both  functions  of  a third  variable,  the  time.  However,  we 
shall  indicate  here  briefly  the  case  of  singularities  of  higher  order.  These 
correspond  to  points  (x1 , y* ) which  are  not  only  common  seros  of  ?(x,y)  and 
<2(x,y),  but  also  stationary  points  for  them.  For  example,  the  equation 

_ 3X2 

P Ty 


admits  the  general  integral 


- x^  ♦ 


This  consists  of  a system  of  cubic s of  which  the  curve  y* 
the  origin. 


x?  has  a 


cusp  at 


The  equation 
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For  an  equation  of  type  (6),  necessary  and  sufficient  condition  of 
integrability  is  that 

P(x,y)  dx  ~ Q(x,y)  dy  - 0 (11) 

be  an  exact  differential,  i.e.. 


IsP  tQ  n 

— ♦ —I  » 0 

) v -r 

v ¥ - — 

The  solution  is  then  expressed  with 

/Xp  dx  ♦ /yQdy  ■ c 
xo  yo 

The  condition  of  integration  is  alvr.ys  verified  when  the  variables  are  sep ar= 
able,  i.e.,  P - P(x),  Q * Q(y),  or  vice  versa. 

If  P(x,y),  Q(x,y)  are  homogeneous  functions  of  x and  y of  the  same 
degree,  the  separation  of  variables  can  be  achieved  with  the  substitution 
y - ax  which  transforms  (10)  into 

[p(1,z)  ♦ z Q(l, z)]  dx  - x Q(l,z)  dz  - 0 

More  generally,  if  (11)  is  not  an  exact  differential,  it  may  be  transformed 
into  one  by  multiplication  with  an  integrating  factor  p.(x,y),  such  that 

n(P  dx  - Q dy)  ■ 0 

is  a total  differential.  There  exists  an  infinity  of  integrating  factors  which 
are  solutions  of  the  partial  differential  equation 

. o 

1y 

i e Q • $ 

P lli  ♦ qIE*|1(J£**£).o  (12) 

"iy  ~bx  "iy  *x 

The  solution  of  (12)  is  in  itself  a problem  more  difficult  than  the  solution 
of  (11),  but  fortunately  ws  only  nssd  to  find  & particular  integral  of  (11) 
and  this,  in  certain  cases,  can  be  found  easily.  As  an  example,  consider  any 
linear  differential  equation  of  first  order 

^ ♦ a.  (x)  y - b.(x)  (13) 

dx  x 

In  this  case  we  have  P « b^  - a^y,  Q - 1,  arid  Eq.  (12)  becomes 


)• 
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(ti 

A * - -icular  -integral  of  this  equation  is  |i  ■ p(x)  obtained  from  = pa, 
i.e.,  <** 

( \ A i(*)d* 

n(x)  - V A 

With  the  use  of  this  integrating  factor,  Eq.  (13)  is  transformed  into  the 
exact  differential 

/a,dx  /a,dx 

e (b1  - a^y)  dx  - e dy  - 0 


whose  general  integral  according  to  (12)  is 


y - y0  - e 


-Ajdx 


JC  f L 

C ♦ / (-a^y  + b)  e dx 

*0 


-/a^  dx  -/a-,  dx  x /a^dx 

y-Ce  * ♦ « ' / b^(x)  c dx 


In  the  case  of  the  nonlinear  Bernoulli  equation 

^ ♦ a1(x)  y ♦ ah(x)  yn  - 0 (Ui) 

dx 

Eq.  (12)  becomes 

(*i  y ♦ • /*)  ^ ^ y”"1)  - 0 

'lx  fix 

A particular  solution  of  this  equation  is  not  easily  found.  However,  Eq.  (lU) 
can  be  transformed  into  a linear  equation  by  means  of  the  substitution  z - yl-R 
and  becomes 

^ ♦ (l-n)(a1z  + an)  - 0 

The  Riccati  equation 

^ y + *2  y2  * bl  (15) 

dx  X A 

which  contains  the  Bernoulli  equation  as  a particular  case,  cannot  be  inte- 
grated by  a general  method.  Its  solution  is  simplified,  however,  if  any 
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(bi  - »,y)  ^ ^ 

1 'ix  1 

A »''\rticular  Integral  of  this  equation  is  (i  - p(x)  obtained  from  2i±  - ^a. 


i.o.. 


dx 


>x( 


/a1(z)dx 

x)  - e 


With  the  use  of  this  integrating  factor,  Sq.  (13)  is  transformed  into  the 
exact  differential 

/a,dx  j\dx 

i (b.^  - a^y)  dx  - * dy  - 0 

whose  general  integral  according  to  (12)  is 

-yv^dx 


x 


y - yQ  - e 


y*a_  dx 

C ♦ / (-a,y  ♦ b)  e dx 

*0  ^ 


i.e.. 


y - C e 


-/a,  dx  -/a,  dx  x /a,  dx 

♦ « / b.(x)  e dx 

X_  ± 


In  the  case  of  the  nonlinear  Bernoulli  equation 

^ ♦ a1(x)  y ♦ a^x)  yn  - 0 (lit) 

dx 

Eq.  (12)  becomes 

(*1  y ♦ *3^)  ^ ^ ♦ ^(»i  ♦ nah  y1"1)  - 0 

'lx  Ox 

A particular  solution  of  tills  equation  is  not  easily  found.  However,  Eq.  (lit) 
can  be  transformed  into  a linear  equation  by  means  of  the  substitution  z - yl“n 
end  becomes 

~ + (l-n)(a1z  ♦ an)  - 0 
The  Riccati  equation 

^ ♦ a,  y ♦ a2  y2  - b1  (15) 

dx 

which  contains  the  Bernoulli  equation  as  a particular  case,  cannot  be  inte- 
grated by  a general  method.  Its  solution  is  simplified,  however,  if  any 
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particular  integral  y*  is  known.  As  a matter  of  fact,  letting  y * y*»  z the 
equation  is  changed  into  the  Bernoulli  type 

— ♦ (2y*a2  + a,)  z ♦ agZ2  - 0 

dx  ■*" 

Passing  now  to  other  possible  forms  of  Eq.  (3),  a solution  can  be 
obtained  at  least  theoretically  when  $(x,y,p)  is  a polynomial  of  degree  n in 
p.  As  a matter  of  fact,  writing  the  equation  in  factorized  form 

*(x,y,p)  “ (p-*1)(p-<2)  ••••  (P— <n)  * °>  (l6) 

one  can  show  that  the  general  integral  is  the  product  of  the  integrals  of  the 
equations  obtained  equating  to  zero  each  factor  of  (16). 

Other  cases  of  interest  are  those  in  which  ♦(x,yJjp)  can  be  written 
in  one  of  the  forms 


X - 

%(p) 

(17) 

y - 

f2(p) 

(18) 

y - 

xf3(p) 

(19) 

The  solutions  of  these  equations  can  be  written  in  closed  parametric  form  in 
terms  of  p.  One  has  respectively 


X - f-^p). 

y - c + / pf'(p)  dp 

(17*) 

y - *2(p)> 

X - C ♦ / - fl(p)  dp 

(18') 

P 

y - xf  (p), 

J r f-»  (p)  dpi 
x - c exp  < / ■?  ILL — - 

> (19') 

3 

[ P - ^3(p) 

Similarly,  when  *(x,y,p)  is  linear  in  x and  in  y,  of  type 

*(x,y,p)  - x^p)  -y  + f(p) 

with  ^ (p)  pp,  the  solution  is  obtained  in  parametric  form  as 

x m fk(p)>  y - ^(p) 

whsrs  f^(p)  is  the  integral  of  the  linear  equation 

dx  _ f'(p) 
dp  p-'fCp) 


x 


P=t(p) 
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and  f^(p)  « xf(p)  +<r(p). 

On  the  other  hand,  the  roots  p^  of 

P - 'f  (p)  “ 0 

satisfy  the  condition  dp/dx  ■ 0 and  furnish  a number  of  singular  solutions  of 
type 

7 - x'f(pi)  + f (p±) 

In  certain  cases  a given  differential  equation  can  be  transformed 
into  an  integrable  form  with  a substitution  of  variable.  For  instance,  the 
transformation  of  Legendre 

x - P,  y - IP  - Y,  p « X 

transforms  $(x,y,p)  » 0 into  the  equivalent  equation 

*(P,  XP-Y,  X)  - 0 


III.  Methods  for  the  Determination  cf  Approximate  Integral  Solutions 

Unfortunately,  in  the  great  generality  of  cases,  a rigorous  solution 
of  a nonlinear  differential  equation  cannot  be  found  and  for  thi3  reason  it  is 
neoessary  to  have  recourse  to  methods  of  approximation.  These  generally  apply 
only  to  a differential  equation  of  type  ( f>) 

p - f(x,y) 

and  in  a domain  D in  which  f(x,y)  is  analytic  in  x and  y.  They  provide  the 
solution  as  a limit  of  a sequence  of  functions,  or  as  3um  of  an  infinite 
series,  or  as  a numerical  expression. 

a)  Method  of  Successive  Approximations  (Iteration) 

If  Xq,  yQ  is  a point  of  the  domain  D,  the  solution  of  (5)  satisfies 
the  integral  equation 

y - y0  + /*  f(^,y(^))  d|  (20) 

where  y(^  ) is  unknown  and  x is  chosen  within  D.  Let  us  consider  the  sequence 

*1  • > y0'< d f 

*0 
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y2  * yQ  + Sx  f(t  .*yi)  d ! 

x„ 


It  can  be  shown  that  this  sequence  possesses  a limit  for  n->oo,  and  that  this 
limit  is  the  desired  solution  of  (5K  Uhile  this  procedure  is  very  general, 
the  convergence  of  the  sequence  obtained  is  rather  slow.  In  order  to  improve 
the  latter,  in  special  cases  modifications  can  be  introduced  in  the  method. 
Suppose,  for  instance,  that  Eq.  (5)  is  of  the  form 

p - fx(x,y)  * f2(x.y)  (21) 

where  e ^2<<^1  wi^hin  ^he  domain  D.  If  a solution  of 

P “ *]_(  *,y) 

can  be  found,  one  assumes  it  as  zero  approximation  and  obtains,  the  first  ap- 
proximation solving 

P - 7 ) ♦ e f2(x,yo) 

Similarly,  the  second  approximation  is  obtained  as  solution  of 

p - f^Xjy)  + e fgtx^) 

and  so  on. 

Another  modification  of  the  method  is  obtained  considering  the  equa- 
tion 


dx  _ 1 

dy  f(x,y) 


which  sometimes  may  be  integrated  more  readily.  One  has  similarly 


x 


W 

% n 


- ±i 

f) 


(22) 


and  correspondingly  the  approximation  sequence  is 


*1 


" *o*/ 

y„ 


f(x_ 


1 
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*2 


d( 

^(3^,  f ) 


b)  Method  of  Integration  in  Series  ( Perturbation) 


While  with  the  previous  method  the  solution  is  expressed  as  limit 
of  a sequence,  it  can  also  be  expressed  as  sum  of  a converging  infinite  series. 
In  a domain  of  analyticity  of  f(x,y)  surrounding  the  initial  pair  (xQ,yQ)  one 
can  compute  all  successive  derivatives  of  y at  (xo^Jo)  by  means  of  the  given 
differential  equation.  Therefore,  by  Taylor  series,  one  has 


y 


(£2a)f  ♦ 

▼ • a • 

2! 


(23) 


4f  1 


r 

I. 


sr 

f. 

5 


w 

1 


r 

k 


I 


V 

i 


In  practice  the  computation  of  successive  derivatives  may  be  rather  cumbersome. 
If  one  accepts  the  principle  tnax  y can  be  expressed  ao  a power  series,  its 
coefficients  can  be  obtained  by  substituting  such  an  expression  with  indeter- 
minate coefficients  for  the  dependent  variable  in  the  differential  equation. 
Since  it  must  make  the  equation  an  identity,  orve.  can  readily  obtain  recurrence 
relations  between  the  coefficients  of  the  powtr  series  by  comparing  terms  in- 
volving the  same  power  of  the  variable. 

The  method  of  solution  in  series  may  be  applied  also  in  the  neighbor- 
hood of  some  types  of  singular  points  of  the  equation,  namely,  those  in  corres- 
pondence of  which  the  solution  is  expressed  in  the  form 

7 - U-x0)r  j»0  + •i(*-*0)  ♦ a2(x'xo)2  * •••} 

with  a07^Q  and  r any  real  number.  To  investigate  the  possible  existence  of 
such  type  of  solution,  one  substitutes  the  series  (2U)  into  the  differential 
equation  and  equates  to  zero  the  coefficient  of  the  terra  of  lowest  degree  in 
(x-x0).  If  this  coefficient  is  independent  of  r,  the  expression  (2h)  cannot 
be  used  to  represent  the  solution  in  the  neighborhood  of  (xo,y0).  Otherwise, 
one  determines  r from  this  equation  and  then  proceeds  to  the  evaluation  of 
all  other  coefficients. 

When  the  differential  Eq.  (5)  can  be  written  in  the  form 

P “ f(*,y)  - fx(x,y)  ♦ « f2(x>y) 

with  t f2«f]_,  it  is  possible  to  express  the  solution  with  an  infinite  series 
whose  terms  are  not  simple  powers  of  the  independent  variable,  but  more  general 
functions  of  such  variable,  simultaneously  defined  in  a certain  interval  of  x. 
Por  this  purpose,  it  is  observed  that  in  the  given  domain  f(x,y)  can  be  con- 
sidered an  analytic  function  not  only  of  the  (complex)  variables  x and  y,  but 
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also  of  e;  therefore,  the  solution  can  be  expressed  in  a ssrias  of  powers  of 
fi,  i » e o , 

y - y0U)  e y1(x)  + e2  y2(x)  + ...  (25) 

Substituting  the  expression  (25)  into  the  differential  equation,  one  obtains 

• n ^ • 

y0  + eyl  + 6 y2  +*“*  “ Vx*Vcyi+--)  + yo+eyl+* ' ^ 

Now,  equating  the  terras  of  same  degree  in  e,  one  obtains  a system  of  recurrence 

equations  which  permits  the  determination  of  the  functions  yQ,  yi,  yj 

For  instance,  yQ  is  the  solution  of  the  equation 

yo  “ fi(x*y) 

This  special  case  of  the  method  of  integration  by  series  is  known  as  "perturba- 
tion method".  A close  analogy  exists  between  perturbation  and  iteration 
methods;  depending  upon  the  particular  differential  equation,  one  or  the  other 
may  he  more  convenient,  but  in  general  both  become  rather  cumbersome  after  the 
second  or  third  approximation. 

In  the  application  of  each  of  the  methods  previously  discussed,  it 
is  necessary  to  satisfy  the  given  initial  or  boundary  condition.  Since  the 
equations  here  considered  are  of  the  first  order,  this  condition  reduces  to 
one  relation,  for  instance,  the  initial  value  y(0).  In  addition,  if  it  is 
known  that  the  equation  possesses  a periodic  solution  (in  response  to  an  ex- 
ternal driving  terra,  or  as  free  solution  of  a degenerate  system),  a periodic 
condition  of  the  type  y(x  ♦ 2n)  - y(x)  must  be  also  satisfied.  A degenerate 
system  is  one  which  is  represented  by  a differential  equation  of  second  order, 
which,  on  account  of  the  relative  smallness  of  the  coefficient  of  the  second 
derivative  or  of  that  of  the  function  y (spring  constant)  reduces  actually  to 
a differential  equation  of  first  order.  A system  of  the  type 

•j2  - f (y)  - 0 (26) 

dx 

possesses  a periodic  solution  if  the  corresponding  path  in  the  plane  (dy/dx,  y) 

ffy)  | 


Fig.  U “ Multivalued  Function 
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is  a closed  curve.  This  implies  that  f(y)  must  be  a multivalued  function 
(Fig.  U)  since,  for  some  values  of  y,  f(y)  will  assume  at  least  two  different 
values  on  the  solution  path.  Thus,  the  function  f(y)  in  general  will  be  non- 
analytical,  including  radicals  or  discontinuities.* 

By  means  of  the  methods  of  approximation  described,  the  solution  is 
expressed  as  an  infinite  series.  In  order  to  satisfy  the  initial  condition 
it  is  sufficient  to  impose  that  the  first  term  of  such  a series  assumes  the 
required  initial  value  and  all  ether  terms  assume  correspondingly  zero  initial 
value.  In  addition,  in  order  to  satisfy  the  periodicity  condition,  when  this 
is  required,  one  must  equate  to  zero  all  the  coefficients  of  terms  in  the  series, 
possessing  a periodicity  other  than  the  required  one. 


c)  Method  of  Sinusoidal  Analysi s 


Wnsn  it  is  known  that  the  solution  presents  periodic  terms  of  a cer- 
tain frequency,  one  can  express  y as  a Fourier  series  in  x and  determine  its 
coefficients  b^  direct  substitution  into  the  differential  equation.  By  this 
procedure,  the  integral  of  the  nonlinear  differential  equation  is  obtained 
from  the  solution  of  a system  of  equations  which  are  not  of  differential  type, 
but  are  of  algebraic  or  transcendental  type  depending  upon  the  nature  of  the 
nonlinearity  existing  in  the  original  differential  equation.  For  instance, 
given 

^ - hy  - f(y,t)  (27) 

dt 

where  it  is  known  that  y is  a periodic  function  of  t,  one  can  assume  in  general 


(a  sin  ncot  + b cos  mot) 
' n n 


y - 


oo 

n-1 


03 

y*  * / (con  a cos  ncot  - on  b sin  ncot) 

~ i n n ' 

n>*l 


There  follows  that 


oo 

f(y>t)  - 2Z  sin  ncot  ♦ f^  cos  ncot) 


w.iere 


n Tf 

~ / f (y,  t)  sin  ncot  dart 
J"r*  n 


/ f ( y r t ) cos  ncot  dcot 

-TX 


* 


Andronov,  Chaikin  - Theory  of  Oscillations,  Princeton  Un.  Press  191*9.  p.  138, 
Chapter  IV. 
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Clearly  the  solution  of  the  equations  obtained  by  identifying  coefficients 
of  isofrequential  terms  is  impossible  in  the  general  case.  If,  however,  it 
is  permitted  to  assume  that  only  a limited  number  of  harmonics  exist b,  then 
the  equations  may  be  simplified  conveniently  to  permit  the  solution  for  the 
values  of  the  various  coefficients  of  the  Fourier  series. 


In  particular  one  might  combine  the  foregoing  procedure  with  a pro- 
cedure of  perturbation*.  For  instance,  one  assumes  that  as  first  approximation 
y can  be  expressed  by  means  of  its  fundamental  term  only,  i.e„ , y * a^  sin  c»t  ♦ 
b^  cos  cot.  Consequently,  Eq.  (27)  reduces  to 

& - hy  ■ f..  sin  cot  1-  f2i  008  «t 
dt 

where 


1 

n 


1 

n 


sin  cot  + b^  cos  cot. 


sin  cot  ♦ b-^  cos  cot. 


*] 

*] 


sin  cot  dcot 


cos  cot  dcot 


Letting  now  y - Re  T,  where  Y • (b^  ♦ ja^)  - A e - A one 


can  write  equivalently 


[*>  - hJ  T - (f21  * ifU>  !>t  ' T 


i.e 


• 9 


(28) 


If  the  nonlinear  term  f(y,t)  contains  also  derivatives,  the  resultant  ex- 
pression N will  be  a function  of  A and  jo.',  i.e.,  N(A,  jco). 


The  term  N(A,  jco)  is  called  describing  function  of  the  nonlinear  ex- 
pression f (y,t) } it  is  seen  that  it  is  a function  of  the  complex  amplitude 
X - Aej*  and  of  the  frequency  co.  If  f(y,t)  is  linear  in  y,  the  corresponding 
describing  function  becomes  the  transfer  function.  In  general,  we  have 

f(y,t)  - Re  N(A,  jco)  Y - p$(y,t)  (29) 

or  identically 

f(y,t)  - j|N(A,  jco)  Y ♦ N(A*  jco)  Y*j  - *i*(y,t) 


* B.  V.  Bulgakov  - Periodic  processes  In  free  pseudo-linear  oscillatory  sys- 
tems, Jour.  Franklin  Inst.,  Vol.  235,  June  19h3,  p.  591-616. 

B.  C.  Johnson  - Sinusoidal  rnalysis  of  feedback  control  systems  containing 
nonlinear  elements.  Trans.  AIEE,  July  1952. 
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The  present  procedure  of  approximate  solution  is  valid  when  in  (29)  one  can 
consider  p*(y,t)  very  small  and  neglect  it  in  first  approximation.  Then  the 
amplitude  AeJ^is  determined  graphically  or  analytically  solving 

J®  - h - N(A,  j«) . 

A better  approximation  to  the  actual  solution  can  be  obtained  when  the  differ- 
ence (29)  is  small  with  respect  to  the  other  terms  in  the  differential  equa- 
tion, by  writing 

^ - hy  - Re  N(A,  jco)  Y - p$(y,t)  (30) 

dt 

In  this  relation  p is  a coefficient  indicative  of  smallness.  Now  following  a 
perturbation  method  we  imagine  to  expand  the  true  solution  y into  a power 
series  in  p,  i.e., 

2 

y - yc  ♦ wri  ♦ y2  ♦ — 

The  first  approximation  is  obtained  as  Indicated  beforehand  letting  p - 0. 

To  obtain  the  second  approximation,  one  differentiates  Eq.  (j0)  with  respect 
to  p and  then  lets  p - 0.  There  follows 

dyi  _ 

— - hyx  - Re  H(A,  Jco)  Y}  - *(yc)  (31) 

dt 

Expanding  *(yo)  in  Fourier  series  one  has 

4-  CD 

$(y  ) - i aJn6st 

o p a_.  n 

c ^OD 

where 

* - - /n*(y  ) dest 

71  -71 

There  follows  for  the  second  member  of  ( 31) 

dy.  ♦cp 

- hy1  - Re  N(A,jco)Y1  -1^1  V.  Jnort, 

dt  2 0 ? -35  n 

n^0,+ 1 


Now  the  condition  of  periodicity  requires  that  $1  - - 0.  This  condition 

is  already  verified,  however,  on  account  of  (29).  Therefore,  the  term  of 
second  approxim&ticn  is  obtained  as 


Re  A^e 


4.-*. 


CO 


n-2 


In* 


jncot 


- 00 


-jncct 
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It  contains  a new  arbitrary  constant  A which  must  be  determined  on  the  basis 
of  the  initial  condition,  as  indicated  previously. 


IV.  Linear  iri.fr’erential  Equations  with  Variable  Coefficients 

The  methods  of  solution  described  for  nonlinear  equations  can  be 
used  also  in  the  case  of  linear  equations.  However,  in  the  latter  case,  since 
superposition  applies,  it  is  possible  to  make  recourse  to  special  methods  of 
analysis,  based  on  the  use  of  the  Green's  function  and  of  linear  transforms. 
For  greater  generality  we  shall  consider  a linear  differential  equation  of 
order  n of  type 

|\n(t)pn  ♦ ...  ♦ »1(t)p  ♦ aQ(t)J  y(t)  -j\m(t)pm  ♦ ...  ♦ bQ(t)~|  u(t)  (32) 

which  we  can  indicate  for  brevity  with 


Mp,t)  y(t)  - k(p,t)  u(t)  - r(t) 


(33) 


where  u(t),  r(t)  are  known  functions  (driving  terms).  In  general,  we 
a set  of  n linearly  independent  relations  in  y( t)  which  represent  the 
conditions.  The  system 

[Up,t)  y - r(t) 


1 


i • 1.  2. . ,n 


also  have 
boundary 


(3U) 


might  not  possess  a solution  not  identically  zero  which  together  with  its  n-1 
derivatives  is  continuous  throughout  the  interval  (a,b)  of  the  independent 
variable.  However,  one  can  always  find  a function  whim  formally  satisfies 
the  system  (3U)  but  violates  at  least  in  part  the  condition  of  continuity.  In 
particular  when  r(t)  - 0,  » 0 (homogeneous  system),  such  a solution  is 

called  Green's  function  G(t,  7r).  Such  a function  is  continuous  and  possesses 
continuous  derivatives  of  orders  up  to  and  including  n-2  when  ait— bj  in  addi- 
tion, it  is  such  that  its  derivative  of  order  n-1  is  discontinuous  at  a point 
t within  (a,b),  and  presents  there  a jump  upward  l/an(,t)|  finally,  it  satisfies 
the  given  system  at  all  points  of  (a,b)  except  at  t - ?. 


If  G(t,U)  is  a solution  of 

Up,t)  y - 0 

u^y)  - 0 i 

then  an  explicit  solution  of  the  nonhomo geneous  system 

L(p,t)  y - r(t) 
u.(y)  - 0 


1. . .n 


(35) 


i - 1. . .n 


(36) 


t?»?< *3| »* fr  -t;  y» :; **,«g.. y - r^^ur./^  ;,.*>}*  | ^ v PV  ■^#<v>»},,'^ip: 
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y{t)  - / G(t,?)  r(r)  d ^ 

a 

The  solution  of  the  more  general  nonhomo geneous  system  (3b)  is 


y(t)  - / G(t/t)r(T)d*P*  1TG|  (*)+...  ♦ JrnGn(+) 

where  G^t)  are  the  UillQuG  solutions  of  the  system 

L(Gi)  - 0 


Uj(Gi)  - 0 

Ui(G1)  - 1 


j - 1 . . i-l,i+l  ..n 


Wien  the  equations  of  the  boundary  conditions  Ui(y)  * are  expressed  in 
diagonal  form 


y^v^(o)  - y„ 


v-0  ..  n-1 


the  latter  procedure  provides  as  functions  Gp(t)  the  solutions  of 

L(p,t)y  - ♦ •••  ♦ ®n«>i+i^0^yoJ  5 “^(O  (UO) 

Therefore,  in  this  case,  the  final  solution  (38)  may  be  considered  as  corres- 
ponding to  a system 

Up,t)y  - r(#  ♦ an(o)yQ$'nal(^)  + Jah(o)y1  ♦ Vl^°^yoJ  i & * 


+ • t » ^ 


[an(°)yn-i + ••• + ai(°>yo]«)(t) 


v±(y)  - 0 


i - 1.  . .n 


In  general,  the  rigorous  solution  of  a variable  linear  system  is  not  known. 
Approximate  expressions  for  it  can  be  obtained  by  application  of  methods  of 
iteration  or  perturbation.  For  instance,  writing  the  system  (3b)  as 

I-1(Pjt)y  - L;,(p,t)y  ♦ r(-fr) 

Ui(y)  " Vi(y)  + *i 
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where  parts  of  a differential  expression  and  of  the  boundary  conditions  have 
been  transferred  to  the  right  hand  side  of  the  equations,  one  can  determine 
uniquely  a sequence  of  functions 

J0(t),  yx(t)  ...  yn(t)  ... 

where  yQ(t)  is  such  that  I^Cy^  is  continuous  in  (a,b) 
and  y^(t),...yn(t)  satisfy  the  recurrence  relations 

Ll(yr)  ■ L2^r.l^  * 

Ui(yr)  - vi(yr_1)  ♦ ri 

One  should  observe  the  relative  arbitrariness  of  choice  of  yo(  t) . Another 
iteration  method*  consists  in  Replacing  the  variable  coefficients  in  L(p,t) 
with  their  mean  values  in  the  range  (a,b)  of  the  independent  variable.,  Eq.  ($3) 
is  then  written  _ 

L(p)y  - k(p,t)u(t)  + AL(p,t)  y 

where 

L(p)  - anpn  ♦ ..  + aQ 

AL(p)  - [V*n(t)]  Pn  * •••  + [V^*!] 

Starting  with  yQ(t)  equal  to  the  solution  of 

L(p).y(t)  - k(p,t)  . u(t) 

one  institutes  recurrence  relations  of  type  (U2)  for  the  Eq.  (1*3). 

The  investigations  of  variable  linear  differential  equations  can 
also  be  carried  out  in  the  complex  frequency  domains.  However,  direct  trans- 
formation of  (32)  leads,  in  general,  to  another  variable  linear  differential 
equation  in  which  order  in  p and  degree  in  s are  interchanged  with  respect  to 
Eq.  (32).  For  this  reason,  in  general,  the  new  differential  equation  is  not 
easier  to  integrate  than  the  original  one.  However,  one  can  find  the  Green's 
function  of  (32)  or  more  generally  the  impulsive  solution  of  (32)  WCtj'fc)  when 
u(t)  « &(t)  and  transform  it  into  the  complex  domain  by  means  of 

H(s,t)  - /tW(t,^)  e“s(t"^dr 


and  vi(yQ)  are  finite, 


i ■ 1 . . .n 


(U3) 


* 


S.  A.  Schelkinoff , M.  C.  Gray  - B.S.T.J.,  Vol.  27,  April  191*8,  p.  350-361*. 
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♦ 

It  may  be  shown  that  th9  general  solution  of  the  systems  (36)  can  be  expressed 
by  means  of 


1 C*J°°  St 

Tit)  - — / H(s,t).U(s)  e ds 

2nj  c~;joo 

where  the  product  H(syt).U(s)  should  be  considered  as  the  frequency  domain 
operation  corresponding  to  the  time  domain  integration 


H(s,t)  - /V,*)  U(*)  d? 
o 

which  is  the  analogous  of  (3?)« 


Recapitulating,  the  procedure  of  solution  of  Eq.  (32)  in  the  time 
domain  consists  of  finding  the  Green's  function  G(t,10  [or  more  generally,  the 
impulsive  response  W(t/t)^  and  then  finding  the  total  solution  y(t)  by  applica- 
tion of  the  superposition  integral  in  one  of  the  forms  (37)  or  (38).  In  the 
frequency  domain,  on  the  contrary,  one  finds  the  corresponding  transform  H(s,t) 
of  the  Green's  function  (or  of  the  impulsive  response)  and  then  obtains  the 
transform  of  the  total  solution  y( t)  by  multiplication  of  the  transforms 
H(s,t)  and  U(s).  In  the  latter  process  one  uses  standard  Laplace  transform 
tables,  considering  t in  H(s,t)  as  a constant  parameter.  It  is  ai so  possible 
to  extend  to  the  analysis  of  variable  systems  many  concepts  familiar  to  fixed 
systems**;  for  instance,  for  thp  case  of  variable  networks  one  can  extend  the 
concepts  of  impedance,  admittance,  gain,  and  various  theorems  of  linear  fixed 
circuit  analysis. 


In  order  to  find  the  function  H(s,t),  which  might  be  considered  as 
the  system  function  of  variable  systems,  one  does  not  have  to  find  first  the 
impulsive  response,  but  can  solve  directly  the  given  differential  equation  in 
the  variable  H(s,t).  For  this  purpose  one  applies  familiar  rule3  of  linear 
differential  operators,  letting  in  (32)  u(t)  - est  and,  correspondingly, 
y(t)  « H(s,t),9s't.  It  is  found 


and  equivalently 


L(p,t)  . H(s,,t)  est  - k(p,t)  eSt 


L(p+s,  t)  . H(s^t)  - k(s,t) 


ihh) 


r 

i_ 


1 Ws,t) 


nl  Un 


P0  ♦ 


r 

1U  s,t) 

p*L(s,t) 

H(s,t)  - k(s,t) 


(U5) 


* L.  A.  Zadeh  - Proc.  Vol.36,  1950,  p.  291. 

**  L.  A.  Zadeh  - J.A.P.,  Vol.  21,  Nov.  1950,  p-  1171. 
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It  should  be  observed  that  Eq»  (U5),  in  general,  is  not  easier  to  integrate 
than  Eq.  (32)  since  one  has  to  use  the  same  methods  of  iteration  or  pertur- 
bation mentioned  previously.  For  instance,  in  applying  the  iteration  pro- 
cedure (ij2)  one  might  choose  as  starting  function  y0(  t)  the  function  Ho(s,t) 
defined  as 

H (a,t)  - 

° L(s,t) 

which  might  be  considered  as  first  approximation  to  H(s,t)  then  the  coef- 
ficients a*(t)  of  L(p,t)  do  not  vary  appreciably  over  the  duration  of  the  im- 
pulsive response.  This  method  corresponds  to  replacing  the  differential  equa- 
tion (U5)  with 

L(s,t)  . H(s,t)  - K(s,t)  - y~  £ P*  H(s,t) 

>1  Jl 

where  all  derivative  terms  have  been  moved  from  the  left  to  the  right  hand 
side.  The  recurrence  Eqs . (U3)  are  correspondingly  written 

n _ . 

L(a,t)  . H(s,t)  - K(s,t)  - > - 6 P*  H .(s.t) 

Jl  ^ 


7.  Some  Examples  of  Application  of  Methods  of  Approximate  Integration 

In  order  to  show  the  technique  of  application  of  methods  of  approxi- 
mate integration  we  shall  work  out  as  examples  the  solution  of  some  differen- 
tial equations.  For  the  purpose  of  checking  the  order  of  approximation  ob- 
tained, we  choose  first  a linear  differential  equation  which  can  be  also 
solved  rigorously. 

Specifically,  we  consider  the  linear  differential  equation 


--  + [a  + jbs(t)J  y - 1 (U6) 

dt 

where  s(t)  - sin  cat,  0<b«ca  . Its  rigorous  solution  is  known 

and  asymptotically  (steady  state)  is  given  by  the  particular  integral 

_ e-at-jb/s(t)ut  j t eat+jb/s(t)dt  dt 

-CD 

Using  the  Fourier  series  expansion 


r J 5 003  (-i)n  J (-)  ejn<ot 


n 'co' 


+ 00 
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one  has 


^ at-1  - cos  cat 


/ e 

-oo 


ca 


^-/t«atZ<-3!n  Vs>  «jMt  « 


- 00  - co 

Therefore,  the  asymptotic  solution  of  (U6)  is  written 
. ♦oo 


3<t>  - ^ 5 cos  “*  ^ <-J)n  •’„<!>  SV3S  <w> 


The  latter  expression  can  be  modified  if  one  applies  the  relation 
b 

cos  m > .m 

■ / 1 c x . 

m co' 


b 

6J  - cos  «t  M y ra  J (bj  ejm*t 


-oo 


There  follows  that  (U7)  can  be  written  equivalently 


y(t)  - Z 1: 


itt-n  . /b>  . /b\  e 


j (jM-n)cot 


.10-11  » ,-UN  /O') 

•'  ‘m'M/  Jn'co' 

m n a ♦ jnca 

This  expression  can  be  ordered  in  the  form  of  an  exponential  Fourier  series 
of  type 


♦co 


y(t) 


-oo 


Qk « 


where 


G,  • — f y(cat)  e~^KW^  dcot 


in  -« 

We  have,  for  instance. 


„ V~  ,n>n  T ,bv  T ,bN  sin(m+n)n 

» ^ ^ ^ " <S>  " <5>  (».).  (a+Jnw) 


ra  n 

Since  in  this  expression 


sin(nH-n)n 

(m+n)n 


if  m+n 


f° 

l-° 


there  follows  that  m ■ -n  and 
♦oo 

o - y~  j 2 

o Z n '•co 

-00 


(U8) 


(U9) 


a+  ,jnco 
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I 


t 

* 


r ; ' 


To  derive  this  formula  we  have  used  the  relations 

j-2n  - (-l)n,  Jn  . j_n  - (-l)n  Jn2 

Equation  (U9)  can  also  be  written 

OD 

2 1 


Q - J 
o o 


1 V , 2 2a 

i *4  " 7TT 

On  the  other  hand,  more  generally  we  have 


TT 

co 


a, 


- -I  Z r Jn(|)  AiSiSSih- 

(a+  jn£o)x(*^n-k) 


2n  -n 


m n 


Since 


sin(m+n-k)n 

(m*n-k)n 

c° 

4 if  m+n-k  n ' 

Cl  o 

there  follows  m+n  - +k 

and 

°k  - 

♦ 00  . 0 . 
X *k'2”Ws> 

-oo 

j rb)  1 

n “ a+jnco 

In  particular 

n 

, JoJl 

00 

A ▼ 

-Wa->)+J: 

Q - j ♦ 3 J 


n 


~1 T~T 

a ♦ co  n 


In  the  coefficients  of  the  series  (U8),  there  appear  Bessel  functions  of  first 
kind,  of  argument  - b/oo  where  b/co  ^<1.  As  a first  approximation  one  can  re- 
place them  with  the  limit  value 


lin  Jn(^) 


sn 

nl  2n 


J (9) 
o • i • 


.1  O'*-* 
n*  >'  ' 


1 

2 


We  shall  now  derive  the  approximate  solutions  of  equation  (U6)  obtained  by 
application  of  methods  of  iteration  and  of  perturbation. 


A 

t 


! 

•h 


1 
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Method  of  Iteration.  We  write  Eq.  (46)  as  follows 

+ ay(t)  ■ 1 - Jb  s(t)  y(t) 
dt 

Operating  the  change  01  variable  bt  - 'i',  there  follows 
d_  d d~C  b d 

dt  <fl?  dt  d-b 

and  Eq.  (50)  becomes 


i.e.. 


b — y(r)  + ay(2*)  - 1 - jb  s(t')  y(-t) 

d r 


— 7 ♦ - 7 - - - J a(T)  y 
dr  b b 


(50) 


By  assumption  | s( t) |—  1 and  we  assume  that  for  our  purposes  the  term  s(t).y 
can  be  considered  small.  Letting 

y(t)  - 70(t)  ♦ yx(t)  + y2(t)  ♦ ...  ♦ yn(t)  ♦ .... 


where  yn(t)  is  of  order  of  |s(t)|n,  substituting  into  Eq.  (50)  and  equating 
terms  of  the  same  order,  wa  have  the  system  of  linear  differential  equations 


• 

v ♦ — y ■ 
-0  b ^0 

1 

b 

(51) 

-i  »(T)  yQ 

(52) 

-i  s(T)  yx 

(53) 

From  Eq.  (5l)  the  steady  state 

solution  is 

v ft3  - i 

-o'  ' 


There  follows  from  (52) 
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i .a.. 


_ i r 


yx(t)  - -3  i-5-/  e H rs('^)  d r 


~ — f — 


Substituting  (5U)  into  Eq.  (53)  > we  have 

. " i 1 

y2  + | y2  - -s(t)l-A_/e  b s(r) 

There  follows 


■ir 


at 


y2(t)  - 


£ ^ A /^*\  A b 


/ .s(?)  / e s(£)  d£  d 


(5U) 


etc.  For  instance,  if  s(t)  - sin  oat  - sin  ^ T,  we  have 

b 


_ “ r 


a >y 

B 


y-, (t)  - -J  / t b sin  ^YdTT  - -j 

b b 


a sin  oat  - o>  cos  oat 

a , 2 2\ 

- (a  ♦ » ) 
b 


There  follows  to  the  first  approximation 


1 . b sin  oat  . to  b 

*tJ  ■ * Ti  - s - 3 -r—y  * 3 — j 


co  b cos  cat 


a ♦ a> 


b(bc  ♦ co  ) 


y2(t)  - 


-b  a b cos  2oot  ♦ 2ocb  sin  2 oat 

.2  . 2\  + «... 2 : 2\ / _2  . r~2T — 


2a(a  ♦ 00  ) 2a(  a ♦ a/)  (a*  ♦ Uaf) 


cob(  ab  sin  2oat  - 2 cab  cos  2 cat) 

+ 77“2 2w  2.2; 

2a(a  ♦ ca  )(a  ♦ Uca  ,! 

Therefore,  the  second  iteration  brings  a contribution  to  the  dc  term  and,  in 
addition,  terms  of  second  harmonic. 


Method  of  Perturbation.  In  order  to  apply  the  perturbation  method  in  Eq.  (U6), 
we  multiply  the  perturbation  term  -jbs(t)  y(t)  by  11  .where  ii  is  a parameter 
indicative  of  the  order  of  magnitude  of  the  perturbation  term.  At  the  end  of 
our  computations  we  shall  replace  it  with  one. 

In  this  analysis  no  study  is  made  of  the  convergence  of  the  series 
obtained  as  a solution.  Following  the  usual  procedure,  it  is  assumed  that, 
even  if  the  series  obtained  are  not  convergent,  the  first  terras  give  a result 
close  enough  to  the  desired  correct  solution. 
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Now  we  expand  y(t)  into  a power  series  in  terms  of  p 

o nl 


t.  : 


Replacing  this  series  into  Eq.  (U6),  we  have 
OL  n . „ 02.  n . oo_ 

y ^ + - yn^  ■ — $»“(*)  yn(b) 

V »i  b Vni  b V —HI 


(55) 


Now  letting  p * 0,  we  have  the  equation  of  the  zero  approximation 

y (TT)  ♦ - y (r)  - i 

Jov  b Jov*-y  b 


(56) 


Differentiating  Eq.  (55)  once  and  then  letting  p - 0,  we  have  the  equation  of 
the  first  approximation 


yx(r)  ♦ - y2(?)  - - J »(*)  y0(^ 

b 


(57) 


Similarly,  differentiating  twice  Eq.  (55)  and  then  letting  p - 0 we  have  the 
differential  equation  of  the  second  approximation 


y2  ♦ ~ 72  m “ 1 2 ®C*)  y^tf) 


(56) 


It  should  be  observed  that  Eqs.  (56),  (57),  and  (58)  are  similar  to  Bqs.  (51), 
(52),  and  (53)  of  the  iteration  method.  The  same  solutions  apply  in  this  case 
also. 


Other  Methods  of  Iteration.  The  given  differential  Eq.  (U6)  may  be  written 
In  the  form 


y(t)  - 


a+3bs(t)  a+jbs(t) 


y (t) 


and  solved  by  iteration  assuming 


y.  - 


n a+ jbs(t)  a+jbs(t)  n’1 


As  solution  y one  assumes 
o 


a+jbs(t) 


t 


T 


s 


&♦  jbs(t)  [a+jbs(t)]3 

For  instance  if  s(t)  « sin  cot,  we  hare 


y2(t)  - 


1b  co  cos  cot 


a+jb  sin  cot  (a+jb  sin  cot)  3 

2 2 2 2 
3co  b cos  cot  - jb(a*jb  sin  cot)co  sin  oat 


(a  ♦ jb  sin  wt; 

Another  procedure  of  approximate  solution  is  obtained  by  writing  the  differen 
tial  equation  in  integral  form 

y(t)  - t - /t(a  ♦ Jbs(t)  y(t)  dt 


Successive  solutions  are  obtained  by  iteration  writing 
yn“t'-/,*0l  + Jbs(t)J  yn-1  (t)  dt 


There  follows 


yQ  - t 


y,  - t - /*  [a  ♦ jbs(t)]  t dt 

x O 

y2  - t - / [a  + jbs(t)]-^t2  - [a  ♦ jbs(t)]  t dt  ^ dt 
For  example,  when  s(t)  - sin  cot  one  has 

t r b 1 

y (t)  - t - / (a  ♦ Jb  sin  cot)  < t-at^i  —x  (sin  cot.  - cot  cos  cot)  > dt 


* 

iZ 
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- t - (a~a2)  ♦ j ^ (1  - cos  cat)  - J ^ (cos  cat  + cat  sin  cot  -1) 


" J (sin  cat  - cat  cos  cat)  + J (sin  cat  - cot  cos  cot)  ♦ 
ca  or 


b . 


fe  ! 

I • 


f ! 

I: 

I 


K 

r 

«! 

I 


2 2 

♦ — (t  - a^~n  2-  — ) - — (sin  2 cat  - 2 cat  cos  2 cat). 

2ca2  2 Uca3 

Another  example  is  provided  by  the  differential  equation 

^ + a(y  ♦ by3)  - A cos  cat 
dt 

which  describes  the  behavior  of  a circuit  consisting  of  a nonlinear  resistance 
in  series  with  a linear  inductance.  If  b^d,  one  can  proceed  by  iteration 
considering  y as  the  limit  of  a succession  and  letting 


dy. 


— - + ay  - A cos  cat  * aby  3 
« n "-1 

If  y(o)  - 0,  one  has 

C • at  A i . ^ \ 

yQ  - — e ♦ ~T~~2  cos  ” w 
a a*>ca 

where  tan  'f  - — , C - - cos  'f 

* a*>ca^ 

Similarly,  y^  is  the  solution  of  the  differential  equation 


dyi 


dt 


+ ayx 


A cos  cat 


ab  4.  - «“at  ♦ cos  (cat  - ?) 

a 


! 


1 


etc. 


w 

m 


A different  procedure  for  the  solution  of  the  given  differential 
equation  consists  in  transforming  it  first  into  the  integral  equation. 


- - sin  cat  - a /*  (y  + by3) 

ca  o 


dt 


\ 


t 


i 

w 

r 


This  equation  can  be  solved  by  iteration  letting 
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y - - sin  cot  - a / (y  + by  T3)  dt 

n co  n n‘-L 


ona  has,  with  the  initial  condition  y(o)  - 0: 


— sin  <ot 
co 

A . ,tl 

— sin  ait  - a J 

oo  n 


— sin  cut  ♦ b(— ) sin^  cot 


00 


CO 


dt 


~ siri  cut  - (1  - cos  cut)  - 


co 


co 


etc* 


~ (-)  < — (1  - cos  cot)  - — ( 1 - cos  3 cat) 

4 co  co  3<o 


VI.  Difference  Equations  and  Methods  of  Numerical  Integration 
Difference  Equations 


Difference  equations  have  particular  importance  for  nonlinear  anaiy— 
3is  because  they  are  used  in  numerical  methods  and  in  some  analytical  pro- 
cedures of  approximate  integration}  in  addition,  they  are  also  used  in  de- 
scribing the  behavior  of  on-off  automatic  control  systems,  tfe  shall  limit 
ourselves  here  only  to  the  treatment  of  lin»er  difference  equations. 

While  in  differential  calculus,  which  deals  with  quantities  varying 
continuously  in  a certain  range,  one  defines  the  differential  operator 
D - d/dt  and  the  successive  operators  D?,  ...  D°,  ..in  the  calculus  of  finite 
differences,  which  deals  with  quantities  varying  dis continuously  in  a certain 
range,  one  defines  the  difference  operator 

* yn  - 3W  - 

find  the  successive  difference  operators 

f*.  ■ ,-•***.  ’ 

In  terms  of  t one  defines  also  the  operator  E - t * 
following  operational  relations 


^ - ? yn 

1,  which  satisfies  the 
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yn  “ 


n+15 


1 


2 

e yr 


- E 7. 


n+l 


- y, 


n+2 1 


etc. 


Line&r  analytical  relationships  among  variables  of  two  corresponding  sequences 
ire  expressed  in  form  of  difference  eqi&tiona  of  type 


QU)  yn  - PU)  vn  (5?) 

wnere  P(a).  Q(A)  are  polynomials  in  A of  type 

QU)  - b0  An  ♦ bj_  A0"1  ♦ ...  ♦ bn 

P(A)  - aft  ♦ a1  Am"x  ♦ ...  ♦ affl 


and  vn  represents  a sequence  of  known  values.  The  general  solution  of  equation 
(59)  can  be  expressed  as  sum  of  the  solution  of  the  homogeneous  equation  ob- 
tained from  (59)  by  letting  P(A)  - 0,  and  a particular  solution  of  (59).  The 
solution  of  a homogeneous  difference  equation  of  type 


(bQ  Lr  ♦ \ * ..  ♦ br)  yn  - 0 

has  the  form 


n 


i-1 


Ai  * 


V 


(60) 


(61) 


where  the  values  of  Tfj  are  the  roots  of  the  transcendental  equation  obtained 
introducing  (61)  into  (60).  One  can  also  express  Eq.  (60)  by  means  of  the 
operator  E and  obtain  a relationship  of  type 


br  y, 


n+r 


r-1  ^n+r-l 


♦ ...  ♦ b.  y„  - 0 
o •'n 


(62) 


The  soi.ution  of  the  latter  equation  has  the  form 


% 


n 

i 


where  ^ are  the  roots  of  the  auxiliary  equation 

b ^ r ♦ ...  ♦ b - 0 

r o 

Because  Eq.  (59)  has  been  assumed  to  be  linear,  the  principle  of  superposition 
applies.  Accordingly,  it  is  possible  to  express  its  solution  by  means  of  an 
integral  expression  similar  to  the  integral  of  Duhamel.  For  this  purpose,  one 
defines  the  unit  step  sequence  Ujj  as 
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Jt. 

jfo 


¥ 


XL. 

as 


lb 

1 

% 

& 

1 

5 

6 

I 

P 


ur  - 0 fear  n - 0 
u 

n 

and  t he  unit  Impulse  sequence  u? 


(63) 


u - 1 for  n »•  0 
n 


t u as 
n n 


u*  >4^*0  for  n i 0 
1 n - 0 


(6U) 


u!  - 
n 


Then  by  application  of  linear  superposition,  one  can  express  any  bounded  func- 
tion yn  in  terms  of  unit  step  or  unit  impulse  sequence  as  follows: 

n 


yn  - y0 


t y u 
Jr  n-r 


or 


Since 


u 

f - "S  y u' 
rn  Jr  n-i 


u 


,-r 


E"r  u* 


_ u , u1  _ 

n-r  n*  n-r  n 

equations  (7)  and  (8)  can  be  written  respectively 


u 

n 


rn  - k * 4 yr  E”’ 

- [y0  ♦ i yx  s"1  ♦ 4 y2  e“2  + .,.J  un  - y(E).Uy 

n 

jr_  - S y E”r  u5  - Ty  ♦ 7i  E_1  u*  - 

n r=U  r n L ° 1 J n 


(65) 


(66) 


(65’) 


— ^ / r\  ..i 

jr  “• 


n 


/ r /■  1 \ 

\.oo  ) 


The  coefficients  of  Y(E)  and  of  Y*(E)  can  be  read  directly  from  a graph  of  the 
sequence  yn  versus  time. 


If  we  now  indicate  with  and  0n;,  respectively,  the  solutions  of 
Eq.  (59)  for  vn  * un  and  vn  - u*n,  the  general  solution  of  the  same  equation, 
can  be  written  as  follows: 


y ** 

Jn 


^S_  „r  1 

L ▼ E A 

r _ n 


i 

1 

i 

\ 


i 

/ 

i 

{ 

\ 


(67) 
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n 


y~  7 E“r  G 
r n 


(68) 


The  functions  An  and  Gn  are  particular  integrals  of  Eq.  (59)  which  can  be 
written  symbolically 

A - Oil  « 


n 


n 


Q(A) 

- u. 

Q(o) 


By  long  division,  letting  A ■ E -1,  and  assuming  that  m<n,  one  has: 

'B  - /«  >•  \®“1 

~ P • • o ~ * -1 

® _ r>  r>  1 t 

'2 


S*2  _ .0  . C,  E-1  ♦ C,  E-2 

Ha)  bo(E-l)n  + b^E-l)"”1  ♦....♦  bn  ° 


Therefore,  An  and  Gn  can  be  evaluated  and  expressed  in  the  forms  (65')  and 
(66'),  respectively.  Prom  these  expressions  one  obtains  the  general  solution 
yn  corresponding  to  the  arbitrary  sequency  vn.  However,  the  latter  can  also 
be  obtained  mors  directly  by  letting  vn  - 7(E) .Un  or  vn  - 7*(E).U'n  in  Bq.  (59). 
One  has  respectively 


P(E-1).V(E) 

Q(E-1) 


u 


_ P(E~1).7*(E)  u. 


n 


Q(E-l) 


Numerical  Integration 

One  of  the  simplest  methods  of  numerical  integration  is  based  on 
direct  use  of  Taylor's  series.  Given 

y*  - f(xsy)  (69) 


with  y(xo)  - y0  one  first  evaluates  from  (69)  y'(*o)*  Now  differentiating 
Eq.  (69)  successively  one  evaluates  y" y' 1 ' „ etc.,  and  correspondingly, 
y'^*©)?  y'  * ' (xq),  etc.  There  follows,  choosing  h small  enough, 

y2 

y(*1)  - y(xQ+h)  - y(xo)+hy'(x0)+  — y' '(»„)♦  •••• 
y'ix-J  - y'(xQ+h)  - y!(xo)-s>hy*»(xo)+  — y'"(x0)+  .... 
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Haring  obtained  y(x^)  and  y*  (x^)  one  evaluates  from  Eq.  (69)  similarly  yM(*i), 
y> ' hxj),  etc.,  and  obtains  y(x2)  ■=.  y(x^  ♦ h).  The  process  is  repeated  until 
the  required  solution  is  obtained.  A check  on  the  accuracy  of  the  computations 
is  made  by  adding  separately  the  terms  of  odd  power  and  those  of  even  power  of 
each  series.  Then,  if  the  series  corresponds  to  x » xn.  the  sum  of  these  two 
results  gives  y(xn+i)  and  their  difference  gives  y^x^^;.  The  latter  value 
should  coincide  with  the  value  previously  calculated, "and  provides  an  indica- 
tion of  the  error  and  a correction  term. 


In  particular  the  method  of  Euler*  as  applied  to  the  differential  Eq. 
(69)  consists  in  using  a Taylor  expansion  limited  to  the  linear  term,  i.e., 


*<*>,>  * • h 

ax  n-1 

An  improved  value  of  (^)  is  found  at  each  st«p  by  multiplying  h by  the 

d*  n-1 

average  of  the  values  of  (dy/dx)  at  the  ends  of  the  interval  x^,  i.e.. 


1 

2 


d*  n-1 


♦ 


A simplified  procedure  of  numerical  integration  is  obtained  replacing 
dy/dx  in  Eq.  (69)  by  an  approximating  polynomial  and  then  integrating  this  over 
any  desired  interval,  Fbr  this  purpose,  one  uses  the  formula  of  Newton 


2 3 

y * - y'  +t  y'  u+  ^ ^ n (u  -»-u)  + ^ ^ n (u^-»-3u^+2u) 
n n /»  / 

c 6 


♦ — —D.  (uii+6u3+llu2+6u) 

2U 


(70) 


where  u -(x-x^h,  and  An  are  difference  operators  of  order  n.  Integrating 

Eq.  (70)  over  the  intervals  x • x ,.x  - x , etc.,  one  obtains  various 

’ ' ' n * n+1'  n-1  n*  ’ 

formulas  for  integrating  ahead  or  for  checking  and  improving  previously  cal- 
culated values.  For  example,  one  has 


r n+1  , . . 

J y'dx  - h 


n 


yl  + 1 ^'y'  + i—  £2y'  + 2 £3y'  * 

J n 2 J n 12  /n  8 * n 720  n. 


f n y'dx  - h 
Xn-1 


y'n-  ± A'y« 


- — A2y'  - — aV’  - — - A*V' 

n 2U 


720 


(71) 


* 


J.  B.  Scarborough,  Numerical  Mathematical  Analysts,  The  John  Hopkins  Press, 
1950,  p.  23U. 
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In  oraar  to  use  these  formulas  one  needs  to  start  the  solution)  since  the 
starting  values  must  be  very  accurate,  this  part  is  usually  very  laborious. 
The  most  common  method  used  to  atari  ihe  solution  is  based  on  the  use  of 
Taylor’s  series,  as  already  described. 


The  approximate  numerical  solution  of  a differential  equation  can 
also  be  obtained  by  replacing  directly  the  equation  with  its  equivalent  dif- 
ference equation.  To  do  this,  one  observes  that  by  application  of  Taylor’s 
series  for  h small  enough 

t h(x) 


y\x+h)  - Tl+hl>*-  D2^...  y(x)  - 

L 21 


and  similarly  in  finite  differences 


There  follows 


y(Xn*h)  - htx^)  - Ey(xn)  - «hDy(xn) 
E - - 1+hD*  — D2+... 


and 


D - i &WS  - ifinil+A)  * l U " £-  ♦ 
h h h [_  2 


(72) 


An  approximate  representation  of  En.  (?2)  in  terms  of  E 


-n 


is 


D 


1 

h 


1 - B~n 

Pn (E_1) 


(73) 


*ho re  ?n(E~^)  is  a polynomial  of  degree  n,  whose  coefficients  might  be  ob- 
tained by  substitution.  For  instance,  it  is  found 


n'-  1 


D 


2 x - E^1 
h 1 + E=1 


n - 2 


n - 3 


D 


r-  8 1-Z  j 

• — ZT T5 r? 

3h  1 + 3E  * 3E  ♦ E 3 


The  errors  of  these  expressions  depend  upon  the  third,  fourth  and  fourth 
difference,  respectively.  Another  less  accurate  representation  of  D which 
is  often  used  is  from  (lii) 


D - - 
h 


(7U) 
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Given  a differential  equation 

Q(D)y  - P(D)v(t) 

one  writes 

y - y(t)  - F(D).T(t)  (75) 

Q(D) 


and  replaces  the  powers  of  D with  the  corresponding  approximate  expressions  in 
E~n.  If  the  highest  order  of  D is  m,  the  approximate  expression  used  for  D 
should  be  accurate  to  the  m-th  difference.** 


One  has 


For  example,  consider  the  differential  equation 
(D  ♦ 1)  y - u(t) 


y - — ■ u(t)  - g—  n un 

1*D  1 ♦ - 1 n 

h 1+E"1 


y « h 
•'n 


1 ♦ E 


-v  u 
“I  n 


2+h°(2«h)E 

This  expression  can  be  evaluated  by  means  of  continuous  division. 


The  application  of  the  previous  methods  to  differential  equations  of 
higher  order  presents  no  difficulty  since  these  can  be  readily  transformed  into 
a system  of  equations  of  first  order. 


The  solution  of  difference  equation  can  also  be  obtained  by  means  of 
the  so-called  relaxation  method.*  The  difference  equation  (given  or  obtained 
by  transformation  from  a differential  equation)  is  written  preferably  in  the 
form  (62),  i0e.s 


b y 
r ■'n+r 


br=l  ^n+r-l 


b y 


- 0 


(&) 
\ ~ / 


This  equation  connects  sets  cf  r consecutive  points.  We  might  assume,  for 
example,  that  y represents  values  of  potential  in  a region  with  given  boun- 
daries. To  begin  with,  one  selects  a rather  wide  square  net  covering  the 
given  region,  and  whose  intersection  points  are  the  values  Xn.  One  fills  the 
entire  area  with  guessed-at  values  of  y,  and  as  a result  the  application  of 
the  Eq.  (62)  now  will  provide  a residual  different  than  zero,  i.e.. 


b y 
r Jn+r 


b , y , ♦ 
r-1  Jn+r-i 


♦ b y 
o Jn 


R(n) 


(76) 


E.  Weber,  Electromagnetic  Fields,  Theory  and  Applications,  J ..  Wiley,  1950, 

p.  260, 

**  Milne-Thomson  - Calculus  of  finite  differences,  London,  1923.  Brown,  B.H.  - 
Math.  Gazette  30  (19U?)»  Iii5  ■ 
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For  the  exact  solution  R(n)  should  be  sere,  and  its  actual  value  depends  upon 
the  initial  value  yn  assumed.  Any  correction  at  n will  affect  the  residuals 
at  all  its  neighbor  points  as  well.  It  requires  therefore  some  little  ex- 
perience to  estimate  the  corrections  needed,  and  it  is  generally  preferable 
to  note  next  to  the  assumed  values  yn  the  residuals  in  brackets.  One  will  use 
the  distribution  of  the  residuals  for  the  second  estimate.  It  is  important  to 
note  that  the  procedure  is  definitely  a convergent  one  even  if  one  starts  from 
a rather  crude  first  guess.  Good  results  are  recorded  by  Strutt*  for  electron 
tube  problems.  The  c-ethod  is  illustrated  in  Cosslett**  and  Zvorkyn,  et  al*** 
for  electron  optical  problems  and  in  Southwell****  for  the  magnetic  flux  dis- 
tribution in  a generator;  many  applications  have  been  made  to  elastic  and  heat 
problems . 


M.  J.  0.  Strutt,  Ann.  d.  Physijc  87,  p.  1$3  (1928).  M.  J.  0.  Strutt, 
Moderns  Hehrgitter-electrcnenrohren,  Vol.  2,  J.  Springer,  Berlin,  1938. 

V.  E.  Cosslett,  Introduction  to  Electron  Optics,  Oxford  University  Press, 
England,  19U6. 

V.  K.  Zworkyn,  et  al,  Electron  Optics  and  the  Electron  Microscope,  J. Wiley 
and  Sons,  19U5>. 

R.  V.  Southwell,  Relaxation  Methods  in  Theoretical  Physics,  Oxford  Univ. 
Press,  19U6. 


